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ABSTRACT 

We report transient quasi-periodic oscillations (QPOs) on minute timescales 
in relativistic, radiative models of the galactic center source Sgr A*. The QPOs 
result from nonaxisymmetric m — 1 structure in the accretion flow excited by 
MHD turbulence. Near-infrared (NIR) and X-ray power spectra show significant 
peaks at frequencies comparable to the orbital frequency at the innermost stable 
circular orbit (ISCO) f Q . The excess power is associated with inward propagating 
magnetic filaments inside the ISCO. The amplitudes of the QPOs are sensitive 
to the electron distribution function. We argue that transient QPOs appear at 
a range of frequencies in the neighborhood of f Q and that the power spectra, 
averaged over long times, likely show a broad bump near f Q rather than distinct, 
narrow QPO features. 
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Introduction 



High frequency quasi-periodic oscillations (QPOs) in the light curves of a ccreting black 
holes (e.g. Morgan et all ( 119971 ): IStrohmayerl ( 120011 ) : iGierlinski et all ( 120081 ): see lRemillard fc McClintock 
( 120061 ) for a review) have attracted attention because they are a potential probe of strong 
gravitational fields. For example, QPO frequencies may be related to / OJ the orbital frequency 
at the innermost stable circular orbit (ISCO) — a feature unique to the strong gravitational 
fields around black holes and possibly neutron stars. Since in general relativity f Q depends 
only on black hole mass M and spin a*, the QPO frequency, together with a mass estimate, 
could be used to infer a*. But this cannot be done with confidence absent a convincing QPO 
model. 

It is difficu lt to discriminate between the many phenomenological and physical models 
of QPOs (e.g. iRemillard fc McClintockl 120061 ) using observat ions alone. Numerical exper 



imen t s offer a potentia l route to test i ng Q PO models (e.g. iKatol 12004 : ISchnittman et al. 



20061 : IChan et al.l 120091 : IHenisev et al.l 120091), but no nu merical models have yet produced 
robust high frequency QPOs. iReynolds fc Milled ( 120091 ). for e xample, find no QPO s at all, 
while in other cases QPOs observed in dynamical variables ( IHenisev et al.l 120091) are not 
obser ved in the emergent radiation calculated from the same simulation (IDexter fc Fragile 
20111). and in still other cases QPOs are observed only at certain times and viewing angles 
dSchnittman et al.l 120061 ). 



Still, accretion disks are complicated radiative and dynamical systems, and faithfully 
modeling all the physics relevant to QPOs is a formidable challenge. In low luminosity 



systems such as Sgr A*, however, the radiative and dynamica' 
and one might hope to build more nearly ab initio models (e.g . 



problems are decoupled 



Moscibrodzka et al. 2009: 



Dexter et al.ll2009l . l2010l : Moscibrodzka et al.ll201ll : iDolence et al.ll2011bh . Searches for QPOs 
in the near- infrared (NIR) and X-ray light curves of Sg r A* itself have resulted in reported 



detections (|Genze 



2006: lMever et al 



et al. 



2003 



Aschenbach et al.l I2004J : iBelanger et al.l 120061 : lEckart et al. 



2006allbl: iTrippe et al.ll2007h . but these may not be statistically significant 



200. 



8|). The observational status of QPOs in Sgr A* is therefore 



(IDo et al.ll2009l : Mever et al 
unclear. 

In this paper we present radiative models of Sgr A* based on three dimensional GRMHD 
simulations. We find clear evidence for QPOs in NIR and X-ray light curves and power 
spectral densities (PSDs) on minute timescales. We tie the QPOs to m — 1 structure 
excited by MHD turbulence in the inner accretion flow. We discuss the QPOs' detectability 
and argue that their amplitudes are likely to be reduced in Sgr A* relative to our simulations 
owing to the likely nonthermal origin of the NIR emission. We argue that the QPOs are 
transient and appear over a range of frequencies near f Q and therefore likely result in a broad 
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bump in PSDs averaged over sufficiently long intervals. Nonetheless, this is the first robust 
detection of a radiative QPO in a self-consistent dynamical black hole accretion model. 



Model 



We model Sgr A* as a hot, optically thin, geometrically thick accretion flow around 
a spinning black hole. The disk orbital and black hole spin angular momenta are aligned. 
We set the dimensio nless spin a* = 1 — _ 2~ 4 ~ 0.94 based on axisymmetric models of the 



quiescent spectrum (IMoscibrodzka et al.ll2009l ). Modeling proceeds in two stages: first we 
evolve the disk using a GRMHD simulation, then we "observe" the disk through relativistic 
radiative transfer calculations. 

The fluid is evolved with the con servative three dimensional GRMHD code harm3d 
dNoble et al.ll2009l : iGammie et al.ll2003l ). The computational volume extends from within 
the horizon to 40 GM/c 2 in radius r, [0.027r,0.987r] in colatitude 0, and [0,27r) in longitude 
cj) with, respectively, 192 x 192 x 128 zones. The zones are regularly spaced in modified 
Kerr-Schild coordinates which are logarithmic in radius and compressed near the equator 
to enhance resolution at small radii and within the main body of the disk. The initial 
conditions consist of a quasi-equilibrium Fishbone-Moncrief torus with pressure maximum 
near 12 GM/c 2 and inner edge at 6 GM/c 2 , small perturbations to the internal energy, and 
a weak purely poloidal magnetic field following isodensity contours. The disk is evolved for 
~ 11, 500 GM/c 3 . After an initial transient phase, at t > 5000 GM/c 3 , the simulation data 
is recorded every 0.5 GM/c 3 for use in radiative transfer calculations. 



The radiation field is evolved with the relativistic Monte Carlo code grmonty ( IDolence et al 



2009 ). grmonty treats synchrotron emission and absorption and Compton scattering. It in- 
cludes all relativistic effects including finite light travel times through the time dependent 
GRMHD simulation data. Dimensional quantities are comp uted from the sca l e free GRMHD 
simula tion by specifying two numbers: M = 4.5 x 10 6 M ( IGhez et al.ll2008l : iGillessen et al. 



2009b. a), and a disk mass unit Ai (effectively the accretion rate). Following IDolence et al. 
( j2011al ). we use a smoothly varying time dependent Ai to remove long term trends in the 
simulation data caused by draining of the initial disk onto the hole. 

The radiation is recorded far from the hole and binned by phot on frequency v in 38 
"cameras" distributed quasi-uniformly over the celestial sphere (see iDolencd l201ll . for a 
detailed discussion of how photons are recorded). The final data set includes broadband 
spectra from radio to 7-rays in each camera wi th an integration time of A = 0.5 GM/c 3 ~ 11 s 
spanning w 5100 GM/c 3 « 31.4 hours (see IDolence et al.l l2011bl . for details). The time- 
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averaged spectra are in general agreement with constraints from sub- mm VLBI, the sub- 
mm spectral slope, and limits on the quiescent X-ray flux (IMoscibrodzka et al.l 120091 ) but 
underproduce NIR flux by a factor of « 10 — 30 (see § 



3. Power Spectra 

The numerical data consists of light curves L"(t, 6, (/)). QPOs are strongest in the plane 
of the disk and absent when the disk is face-oiJE For simplicity, then, we consider only light 
curves at 9 = 7r/2, (j) = k{2n /N c ) ) < k < N C) recorded by the N c = 12 equidistant cameras 
in the equatorial plane. QPOs are detectable at all wavelengths < 100 /xm, but we restrict 
attention to the NIR (3.8 fim and 1.7 fim) and X-ray (integrated from 2 — 8 keV) light curves 
since these are of greatest interest for Sgr A*. 



Following iPress et al.l ( 119921 ). the light curves are divided into A^iocks = 19 equal blocks. 
We compute power spectral densities (PSDs) for each of the 7V seg = A^iocks — 1 segments 
composed of consecutive blocks and average these to obtain smoothed PSD estimates. We 
compute PSDs in both azimuthal wavenumber m and in temporal frequency f n = n/N t A. 
Before Fourier transforming, the data is mean subtracted and Hamming windowed in time. 
The discrete Fourier transform is 

TvwTvw ( _ N j2< m <N c /2 

C=EE w * L u ~ (L"))e-^W N °-<*' N * J /o ~ ^ - at C ; o (1) 

k=0 1=0 I 

and the normalized PSD is 



N t /2 <n< N t /2 



P{v,m,n) = ^^U mn LZ n (2) 



where W ss = N t ^ « o wf is the window function squared and summed as in |Press et al.| 
( 119921 ) and < Ll > is the mean squared signal. Since L v kl is real we need only consider the 



one-sided PSD defined at f n > 0. With this choice m > (m < 0) components circulate 
in the (—(/)) direction. Summing over m is equivalent to summing the PSDs from each 
camera, which is justified by symmetry. Averaging the PSD estimates over segments and 
cameras improves the signal-to-noise-ratio. 

The normalized power spectra are shown in Fig. [H The total power is shown as a heavy 
solid line in the 3.8 /im, 1.7 /xm, and 2 — 8keV bands, and the contributions from each m are 



1 See IShiokawa et al.l ()2Q1 lh for a discussion of convergence of our GRMHD and radiative simulations. 
2 QPOs are detectable simultaneously at all viewing orientations except face-on. 
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shown as lighter colored lines. The ISCO orbital frequency f is shown as a vertical dotted 
line. The main result of this work is the existence of several peaks superposed on a broad 
bump in the power spectrum near f Q in all three bands. We identify these peaks as QPOs. 



The PSDs sho w a power law dependence P 



Mever et al. 


2008; 


Do et al. 


20091)1 



" 2 at low frequency, consistent with 
This low frequency component is 
dominated by fluctuations in the m = (axisymmetric) component of the flow. Near the 
ISCO frequency the m = 1 component dominates the fluctuations and produces a broad, 
Q ~ 1 (Q = v/ Avfwhm), bump in the PSDs. Superposed on this broad feature are three 
peaks with centroids fi = 0.106 min -1 , f 2 = 0.141 min -1 , and f% = 0.166 min _1 @, which may 
be compared to f Q = 0.112 min -1 . The peaks have Q ~ 10. At higher frequency P ~ f~ 3 
and each m > 1 contributes a broad peak near / ~ mf . At the highest frequencies the 
PSDs are dominated by white noise due to Monte Carlo shot noise. The power in m < 
components is negligible. 



4. Discussion 



Where do the QPOs originate? The NIR and X-ray flux is dominated by emission 
from 1.5 < rc 2 /GM < 2.5 close to the equatorial plane. The QPO is therefore generated 
close to the ISCO, at rjsco ~ 2.044 GM/c 2 ) and not far from the event horizon at r^ar ~ 
1.348 GM/c 2 . It clearly probes the strong gravitational field regime. 

The NIR flux is — in our model — thermal synchrotron emission from electrons in the 
high energy tail of the distribution function. The NIR emissivity is sensitive to magnetic 
field strength and temperature, which drop sharply with increasing radius. The mean NIR 
emission is therefore confined to a ring bounded at large radius by declining emissivity and 
at small radius by gravitational redshift and photon capture by the black hol^E 

The mean X-ray flux is dominated by synchrotron photons from the high energy side 
of the synchrotron peak that are Compton upscattered once in almost the same ring that 
generates the NIR emission. This ring is bounded at large radius by a decreasing probability 
of a large energy amplification scattering due to decreasing temperature (0 e = kT e / (m e c 2 ) ~ 



3 Our power sp ectrum does not sample frequencies low enough to see the PSD break reported by 
Mever et all (|2QQ9h . 



4 / 3 ~ (3/2) /i, but the significance of this is unclear 

5 Optically thick disks, by contrast, generate emission in a broad annulus at « 2risco', perhaps this 
explains the lack of high frequency QPOs in the high, soft state in black hole binaries. 
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r -3 / 2 ) and declining optical depth, and like the NIR emission is bounded at small radius by 
gravitational redshift and black hole photon capture. Compton scattering occurs over a 
finite radial range, however, and one would expect that finite light travel time effects would 
tend to average away fluctuations on timescales less than the light-crossing time. This is 
consistent with Fig. [Q where the power at / > f Q drops off more sharply in the X-ray band 
than in the NIR: Compton scattering low-pass-filters the NIR signal. 

The fluctuating component of the NIR and X-ray flux need not arise in the same place 
as the mean signal. We investigated the origin of the variable component by masking out 
emission from r < risco and recomputing segments of the light curve with a strong QPO 
signal. In these segments the QPO disappeared, confirming the importance of sub-ISCO 
emission in generating the QPO. 

4.1. Camera footprint 

Because NIR and X-ray photons are generated so close to the event horizon, gravitational 
lensing, Doppler beaming, and time delays play a role in generating the observed signal. To 
investigate the geometry of the emission region we have focused on a small segment of the 
light curve where a QPO is particularly strong. Part of this curve is shown in Fig. [2 Within 
the segment we generated a vertically integrated map of the radius and azimuth (relative to 
the observing camera) of the origin of photons that are detected in the NIR. The resulting 
map of dN/dxdy is shown in Fig. [31 The disk orbits counterclockwise around the hole and 
the camera is at the far right. Evidently each camera detects emission from a region that 
forms a slender, leading spiral around the hole. The spiral footprint extends for almost 3tt 
radians. 

The camera footprint is radially narrow but azimuthally extended. Since the obser- 
vation process convolves disk structure with the camera footprint, the footprint filters out 
fluctuations from azimuthally narrow structures (evident in Fig. @J indeed most studies show 
flat spatial power spectra for disk turbulence from m = 1 up to m ~ R/ Hx a few). The 
decline in power above ~ 2f Q (and thus, according to Fig. [U with increasing m) is partially 
due to this smoothing effect of the camera footprint. The radial narrowness of the footprint 
also implies sensitivity to radially narrow structures. Radial infall of an axisymmetric emit- 
ting region through the footprint would yield variability. No nonaxisymmetric structure is 
required. Nevertheless, our analysis shows that the QPOs stem from m = 1 structure in the 
emitting plasma. 
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4.2. Underlying flow structure 

What is the m = 1 structure that generates the QPOs? We investigate this by suppress- 
ing nonaxisymmetric structure in each fluid variable in turn and recalculating a segment 
of the light curve with a strong QPO signal. This procedure reveals that the QPOs are 
mainly generated by variations in magnetic field strength Fluctuations in the temperature 
contribute, but at a lower level. Lightcurves in which nonaxisymmetric structure in the 
magnetic field and temperature are removed show little variability. 

If we ignore the dips in the power spectra between the peaks at fi and / 3j the QPO 
is a broad feature with Q ~ 1. The individual peaks at / 1? / 2 , and / 3 have Q ~ 10. On 
average they are, therefore, not very long-lived or coherent structures, as one might expect 
if they are due to turbulence in the disk. They are not orbiting blobs; we see no evidence 
for coherent, orbiting "hot spots" in the disk. 

Figure @] shows a single snapshot of pb 2 in the midplane close to the black hole, which 
shows that there is large-scale structure in the magnetic field. Visual inspection of animations 
of similar images shows that trailing spiral magnetic filaments propagate inward and move 
approximately with the fluid velocity. That the QPO-generating features track the fluid 
velocity strongly suggests that the QPO frequency should scale with the orbital frequency 
at the ISCO, and therefore vary with a*, but we have not yet analyzed models with different 
spin. 



4.3. Significance and persistence of individual peaks 

Are our QPOs statistically significant? To check, we compute PSDs of the full unseg- 
mented light curves and fit each with a sequence of nested models including from zero to 
three QPOs. Likelihood-ratio tests then indicate that all three QPOs are highly significant 
in all three power spectra, with p- values < 5 x 10 -8 . The fitted PSDs and the best-fit models 
are shown in Fig. [51 

Though significant, the QPOs are not persistent. Figure [5] shows the PSDs from the 
first and last half of the simulation, averaged over segments as in Fig. [TJ We see that the 
PSDs from the full simulation are very similar to those from the last half alone. This is 
because the last half of the simulation includes weak flares that dominate the total power. 



6 Using an angle- averaged emissivity barely alters the QPO amplitude, so variations in field direction are 
not essential. 
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We also see that the PSDs from the first half show QPO-like features but with slightly 
different frequencies, though close to f Q . This suggests that PSDs computed from longer 
time series would have broad Q ~ 1 bumps near / OJ but finite length realizations of the light 
curves are likely to show narrower QPO signals at one or more frequencies that may reflect 
real quasi-periodic behavior. Our model shows transient QPOs with Q ~ 10 that appear 
stochastically at frequencies near f Q . 



4.4. Detectability 

Are the QPOs detectable? Recall that our mod el underproduces the NIR flux of Sgr A 



by a factor of 10-30. Also, the observed vF v oc z/ 5 (IHornstein et al.ll2007l ; iTrap et al.ll201lh . 
whereas our model predicts a red spectral slope. This can be understood if the electron 
distribution function has a quasi-thermal core and a nonthermal, high-energy tail. For a 
power-law distribution of electrons with number density n e , dn e /d^ ~ r )~ p ) the emissivity 
is oc z/ - (p -1 )/ 2 5 so the observed NIR slope implies p ~ 2, a common slope for synchrotron- 
emitting sources. The change in NIR emissivity may affect the detectability of the QPO. 

We hav e not yet in corporated a nonthermal tail of electrons into our time-dependent 



model (see iLeungl [2010, where time-averaged models with a nonthermal tail can produce 
Sgr A*'s NIR flux). Naively, one would expect a nonthermal model to be less sensitive to 
B than a thermal model: a nonthermal component with fixed n e has j u oc 5^ +1 ^ 2 , i.e. 
dhijjs/dhiB ~ 1.5, whereas in the NIR our thermal model has j v oc exp(— (V/za^) 1 / 3 ), where 
^ s (O e , B) is a characteristic frequency for synchrotron emission. For parameters appropriate 
to the NIR in Sgr A*, din j u /d\n B ~ 3, suggesting a reduced amplitude for a nonthermally 
generated QPO. On the other hand, the density of nonthermal electrons may be sensitive to 
B and e . That some nonlinear sensitivity of the emissivity to £>, O e , or n e is required is 
consistent with the absence of a QPO near the sub-mm peak, where the emissivity is only 
weakly sensitive to all three. An accurate assessment of QPO strength for any particular 
nonthermal component will require a full, time-dependent model. 

Nonetheless, a prominent feature in the power spectrum close to f Q in geometrically 
thick, optically thin accretion flows seems robust. Hot flows concentrate their emission in a 
narrow ring near the ISCO whereas optically thick disks do not. The emission is (strongly 
or weakly) variable due to the inevitable presence of turbulence in the accretion flow, and 
variation at frequencies / 3> f due to small scale structure will inevitably be averaged away 
by the observing process. 
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5. Summary 

We have performed GRMHD and radiative simulations of the accretion flow in Sgr 
A* and found prominent features near the ISCO orbital frequency in both NIR and X-ray 
light curves. These features are not present for face-on observers or in the millimeter or 
submillimeter. We have shown that: (1) the features have an m = 1 structure on the 
celestial sphere as seen from the source, and are therefore due to m — 1 structure in the 
source; therefore (2) full 2n azimuthal domain models are required to accurately model the 
light curves of similar sources; (3) the variable emission arises near and inside the ISCO and 
therefore probes a strongly relativistic regime close to the event horizon; (4) observations in 
the NIR and X-ray bands are sensitive to a narrow spiral footprint on the disk midplane; 
(5) the variability for our emission model, is dominated by variations in the magnetic field 
strength, with a lesser contribution from variations in the disk temperature; (6) the varying 
features move approximately with the fluid velocity, therefore (7) the centroid frequencies 
should be sensitive to the ISCO orbital frequency and therefore black hole spin. Discovery 
of similar features in the variability spectrum of Sgr A* would be an exciting opportunity to 
probe the spin of the galaxy's central black hole. 
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their hospitality. 
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Fig. 1. — Power spectra of light curves at 3.8 /im (top), 1.7/xm (middle), and integrated 
from 2-8 keV (bottom). QPOs are clearly seen near f Q in all three power spectra, where the 
power is dominated by m = 1 structure. The dotted vertical line shows the ISCO frequency. 



4 



3.8 ]im - 




Q I l I l I l I l I l I l I l I 

8230 8240 8250 8260 8270 8280 8290 8300 

t [GM/c 3 ] 



Fig. 2. — NIR and X-ray light curves over a selected interval showing clear quasi-periodic 
structure. Times shown correspond to the times of detection, ~ 100 GM/c 3 after emission. 
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Fig. 3. — Time- averaged distribution function dN/dxdy for the origin of NIR photons de- 
tected by an observer far away along the +x-axis. The black hole spin and disk orbital 
motion are counterclockwise and the ISCO is indicated by the solid white line. 
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Fig. 4. — A snapshot at t = 8120 GM/c 3 showing pB 2 in the equatorial plane clearly 
indicating the presence of m — 1 structure. 
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Fig. 5. — The noisy central line in each panel shows the PSDs from the full light curves 
without averaging over segments. The smooth lines show the best fit models. The bottom 
and top lines in each panel show the PSDs from the first and last half of the simulation, 
respectively, averaged over segments as in Fig. [U renormalized to account for the lower 
frequency resolution, and then offset by factors of 30. 



